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We consider the statistics of light amplitude fluctuations for the propagation of a laser beam 
subjected to multiple filamentation in an amplified Kerr media, with both linear and nonlinear 
dissipation. Dissipation arrests the catastrophic collapse of filaments, causing their disintegration 
into almost linear waves. These waves form a nearly- Gaussian random field which seeds new fila- 
ments. For small amplitudes the probability density function (PDF) of light amplitude is close to 
Gaussian, while for large amplitudes the PDF has a long power-like tail which corresponds to strong 
non-Gaussian fluctuations, i.e. intermittency of strong optical turbulence. This tail is determined 
by the universal form of near singular filaments and the PDF for the maximum amplitudes of the 
filaments. 

OCIS codes: (190.0190) Nonlinear optics; (190.5940); Self-action effects; (260.5950) Self-focusing 



The self-focusing and multiple filamentation of an in- 
tense laser beam propagating through a Kerr media has 
been the subject of intense research since the advent of 
lasers [J [2j. Multiple filamentation has many applica- 
tions ranging from laser fusion [3|-[5j to the propagation 
of ultrashort pulses in the atmosphere [g, lZJ- Here we 
study the statistics of multiple filamentation, which can 
be viewed as an example of strong optical turbulence with 
intermittency [8J], i.e. strong non-Gaussian fluctuations 
of the amplitude of the laser field [9| . 

Long non-Gaussian tails of the PDF of light amplitude 
fluctuations have been previously observed in filamenta- 
tion experiments [l0| and optical rogue waves [llj . Long 
tails were obtained in solutions of the complex Ginzburg- 
Landau equation with a quintic nonlinearity [12j . The 
analytical form of a long tail of PDF for velocity gradi- 
ent dominated by near-singular shocks was obtained in 
solutions of the forced Burgers equation |l3( . 

Here we describe the propagation of a laser beam 
through the amplified Kerr media by the regularized 
nonlinear Schrodinger equation (RNLS) in dimensionless 
form, 

id z ip + (1 - iae)V 2 V> + (1 + ice)\^\ 2 ip = iebip, (1) 

where the beam is directed along z-axis, r = (a;, y) are 
the transverse coordinates, ip(r, z) is the envelope of the 

electric field, and V = ( -§-, -M- J . The term with a de- 
scribes a wavenumber-dependent linear absorption, the 
term with c corresponds to two-photon absorption, and 
b is the linear gain coefficient. It is assumed that the 
wavenumber-independent part of linear absorption is in- 
cluded in b so that the a-term corresponds to the expan- 
sion of the general absorption coefficient near the carrier 
wavenumber of the laser beam in the Fourier domain. 
Here, e<l and we generally assume a, c, b ~ 1. RNLS 
([T]) is also called the complex Ginzburg-Landau equation. 
RNLS ([1]) can be realized experimentally in numerous 
systems, including e.g. the propagation of a laser beam 
in a ring cavity with a thin slab of Kerr media and ampli- 
fication. In this case the nonlinear phase shift of the laser 
beam at each round trip is small so we can obtain ([T]) in a 
mean-field approximation [14l Il5| with z corresponding 



to the number of round trips in the cavity. RNLS (p} 
also describes the multiple filamentation of an intense, 
ultrashort laser beam in a Kerr media if we average over 
the temporal extent of the pulse [iSj. E.g., multiple fil- 
amentation experiment [171 ] with two-photon-dominated 
absorption corresponds to ce ~ 0.025. 

Neglecting dissipation and amplification, we recover 
the nonlinear Schrodinger equation (NLS) 



id z ip + V 2 Tp+ |^| 2 V> = 0. 



(2) 



NLS describes a catastrophic collapse (also called wave 
collapse) of filaments, max)?/' = \tp\max — ► oo, in a finite 

r 

distance along z, if the optical power TV = f \ijj\ 2 dr is 
above the critical power N c ~ 11.701 [2j. 

The optical power is not conserved in RNLS ([l} for 
e ^ 0. If b > 0, the amplification term on the right 
hand side of ([l| results in an increase of N. If b = we 
assume that N > iV c for z = (e.g. N/N c ~ 10 4 in 
[l2] ) ■ In both cases the modulational instability [2J leads 
to the growth of perturbations of the beam and seeds 
multiple collapsing filaments. These two cases are called 
forced and decaying turbulence, respectively, in reference 
to turbulence in the Navier-Stokes equations Q. The 
statistical properties of these two cases are similar, pro- 
vided in the decaying case we consider distances along z 
at which a relative cumulative decay of N is small. 

In this Letter we focus on a forced case in which a 
dynamic balance is achieved between the pumping of op- 
tical power into the laser beam and dissipation. Figure 
[lja) shows the evolution of N(z) obtained from a numer- 
ical solution of RNLS ([1]). The optical power grows until 
reaching a statistical steady-state corresponding to fully 
developed optical turbulence with N ~ 1200. In this 
regime the amplitude \ip\'is characterized by the random 
distribution of filaments in r and z, as seen in the snap- 
shot of \ip\ for a fixed z in Figure [2] Dissipation is impor- 
tant only when the amplitude of each collapsing filament 
is near to its maximum (see Figure [ljb)) and as well as 
for large wavenumbers k. When \ip\ ma ,x(z) goes through 
a maximum, N experiences a fast decay due to dissipa- 
tion. The influence of periodic boundary conditions on 
the statistical properties of optical turbulence can be ne- 
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FIG. 1: (Color online) (a) N(z) from simulation of RNLS © 
with b = 20, o = c = 1, e = 0.01. (b) The zoom of N(z) in 
a smaller interval in z (solid curve, scale on the right) super- 
imposed the |^| m oa!(^) (dashed curve, scale on the left). All 
simulations used the fourth-order pseudo-spectral split-step 
algorithm on —12.8 < x,y < 12.8 with periodic boundary 
conditions at resolution 4096 x 4096 grid points. Initial condi- 
tions were a superposition of 100 randomly placed real-valued 
Gaussians with amplitudes and radii on [—2, 2] and [1, 2]. 




FIG. 2: (Color online) Snapshot of \ip\ (vertical axis) vs. spa- 
tial coordinates (x,y) for simulation of Figure [l] 



glected if the simulation domain is large enough, so that 
N > N c (in Figure [U N ~ 1007V C ). 

The evolution of each collapsing filament is well 
approximated for lar ge \ip\ by a self-similar radially- 
symmetric solution 0, LL8[ : 



Mt,z)\ = j^Ro(p), P= Wy 



r = r 



(3) 



where L(z) is the transverse spatial scale of a filament 
and Rq is the ground state soliton solution of NLS (0), 
given by \7 2 Ro — Ro + Rq = 0, and corresponding to the 
critical power, N c — J* i? (r) 2 dr 0- 

If e = then L{z) ~ (2n) 1 / 2 (z ~ z) 1 / 2 /{\n\\n(z - 
z)!) 1 / 2 describes a singularity (catastrophic collapse of a 
filament) as z — > z 0, 13 ■ For e ^ the collapse is reg- 
ularized and IV'Imazfz) achieves a maximum \i\)\ r 
at some z — z max . A function 



\maxmax 



dL 
dz ' 



(4) 



changes slowly with z compared to L at z < z max . In the 
vicinity of a collapse, the forcing term in the right hand 



FIG. 3: (Color online) Dependence of |^/)|max(^) for multiple 
individual collapsing filaments (a) in the non-rescaled units; 
(b) the rescaled units (see the text for the description of rescal- 
ing). Individual filaments are extracted from simulation of 
RNLS HI with parameters of Figure [T] 



side of RNLS ([T]) can be neglected; the resulting equation 
can be written in rescaled units z\^\ 2 naxmaxl v\^\ maxmaxi 
and i()/\ip\ m axmax- (Here, we have also shifted z max to 
z = 0.) As shown in Figure EH \ifj\ ma x{z) rescaled in 
these units exhibits a universal behavior for all near- 
singular filaments, — even for e ^ 0, and independent 
of the complicated structure of optical turbulence. This 
universality is a characteristic feature of two-photon ab- 
sorption term in RNLS ([T]), but may not hold for other 
types of absorption. 

Once the amplitude of a filament reaches its maximunr, 
the amplitude decreases and subsequently the filament 
decays into outgoing cylindrical waves as seen in Figure 
O Superposition of these almost linear waves forms a 
nearly random Gaussian field and seeds new filaments. 
Figure 2] shows the probability V(h) for the amplitude 
|^| to have a value h, determined from simulations as 



V(h) = 



J ^(^(r, z)\ — h)drdz 
J drdz 



(5) 



Here, the integrals are taken over all values of r and all 
values of z after the turbulence has reached the statisti- 
cally steady state. We observe that the fit to the Gaus- 
sian distribution works very well for \ip\ < 2 which corre- 
sponds to almost linear waves, while for \ip\ > 3 the PDF 
has a power law-like dependence indicating intermittency 
19]. 

We now show that the power-like tail of V(h) re- 
sults from the near-singular filaments. This approach 
dates back to the idea of describing strong turbulence in 
the Navier-Stokes equations through singularities of the 
Euler equations |9j. Unfortunately, this hydrodynamic 
problem remains unsolved. The forced Burgers equation 
remains the only example of an analytical description of 
strong turbulence in which the tail of the PDF for nega- 
tive gradients follows a well established (—7/2) power law 
13], dominated by the dynamics of near-singular shocks. 
First we calculate the contribution to the PDF from 
individual collapsing filaments. As shown in Figure |31 
the filament amplitude l^lmaa: reaches the maximum 
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FIG. 4: (Color online) (a) V(h) for \i/j\ = h (solid curve) 
for the same simulation as in Figure [TJ Dashed line shows 
fit to the Gaussian distribution and dotted line shows h~ 6 
power law. Circles correspond to the solution of (0. 
(b) H max {hmax) for \ip\maxmax ■ 2231 collapse events with 
IV'Imarrmaa; > 10 are included in simulations in log-log scale 
(red dots). Dotted line shows h~ x power law. Scattering of 
data points for h > 110 is due to lack of statistical ensemble 
for large collapses and is reduced for larger simulation times. 



\ip\maxmax = h max at z = z max , and rapidly decays for 
z > z max . While neglecting the contribution to V(h) 



from z > z n 



we calculate the contribution of an indi- 



vidual filament to V(h) through the conditional proba- 
bility V (h\h ma x) using (O, (|5|), and ([J| as follows 



V{h\h max ) oc 



dz I drS h — 



1 



-Rq 



L(z ma<c ) 



oc / dp p 
dpp 



dLL 3 
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L(z) 



L(z) »\L(z) 



Mp) 



(6) 



Ro(0) 



(~f)h 5 " uvr ' J " \L(z max ) 
= Const h~ 5 6 (hmax - h) , 

where h max = Ro(0)/L(z max ), and Q(x) is the Heaviside 



step function. Here, we have changed the integration 
variable from z to L and approximated 7(2) under the 
integral by its average value (7) as 7(2) — (7) ~ —0.5. 
This approximation is valid for z < z max outside the 
neighborhood of z = z max . 

As a second step we calculate V(h) by integration over 
all values of h max using equation ([6]) as follows 



' \">) / U'lmax ' \P j \'^rnax)irnax\P j rnax) 

= Const h 5 H max (h), (7) 

where V ma x{h m ax) is the PDF for h max = \ip\ ma xmax 
and H max (h) = f^° V ma x( h max )dh max is the cumulative 
probability that \ip\ ma xmax > h. 

Figure EJb) shows H max (h max ). Circles in Figure HJa) 
show the prediction of equation ^\ with H max (h max ) 
from Figure 0Jb). The constant in equation was cho- 
sen to fit the circles and the solid curve in Figure HJa). 
The very good agreement between these two curves, first, 
justifies the assumptions used in derivation of the equa- 
tion ([7]). and second, shows that the intermittency of 
optical turbulence of RNLS ([J) is due to collapse dy- 
namics, which is the main result of this Letter. Figure 
[4] also shows that Hmax{h m ax) is not well approxinrated 
by h m \ x (and is not universal because it depends on the 
parameters a, b, c and e.) Consequently, h~ 6 is only a 
crude approximation for V(h). 
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